function gS = inverseGradient(X)

D = size(X, 1);
invX = inv(X);

gS = zeros(D, D, D, D);
for i = 1:D
    for j = 1:D
        gS(:, :, i, j) = -invX(:, i)*invX(j, :);
    end
end
